Chaotic and Regular Motion in Dissipative Gravitational Billiards 
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We consider the motion of a particle subjected to the constant gravitational field and scattered 
inelasticaly by hard boundaries which possess the shape of parabola, wedge, and hyperbola. The 
billiard itself performs oscillations. The linear dependence of the restitution coefficient on the particle 
velocity is assumed. We demonstrate that this dynamical system can be either regular or chaotic, 
which depends on the billiard shape and the oscillation frequency. The trajectory calculations 
are compared with the experimental data; a good agreement has been achieved. Moreover, the 
properties of the system has been studied by means of the Lyapunov exponents and the Kaplan- 
Yorke dimension. Chaotic and nonuniform patterns visible in the experimental data are interpreted 
as a result of large embedding dimension. 
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From the mathematical point of view, billiards consti- 
tute an interesting class of dynamical systems because 
they exhibit - despite their simplicity - a variety of non- 
linear phenomena, including both regular tori and com- 
pletely chaotic, dense trajectories. Some of them are 
quite realistic and have a direct physical importance. An 
example of such a system is the gravitational billiard in 
which a point-mass particle bounces within a container 
of a given shape and its motion between the bounces is 
not free but ballistic. Obviously, the dynamics depends 
on the billiard shape. For the parabolic boundary the 
system is integrable and all orbits are regular and stable. 
Historically, the first study on the gravitational billiards 
has been performed for the two-dimensional wedge, de- 
fined as two intersecting straight lines 0, Q ■ The motion 
in the wedge is fully chaotic if its vertex angle is larger 
than 7r/2 [3|, otherwise there is a coexistence of regular 
and chaotic behaviour. Those findings were successfully 
tested experimentally The chaotic, as well as regular, 
dynamics is present also for the hyperbolic shape of the 
gravitational billiard which involves both the wedge and 
the parabola as its asymptotic limits .5]. The chaotic 
component prevails for shapes close to the wedge. 

All the above approaches assumed that the collisions 
between the particle and the billiard boundary were elas- 
tic. It is natural to require that, in order to make the 
problem more realistic, one should take into account the 
energy loss and allow for the energy exchange between 
the particle and the wall. Such handling of the dissipa- 
tion is known in nuclear physics as the wall formula [6j, 
derived in the framework of the liquid drop model. The 
atomic nucleus can then be modeled as a billiard pos- 
sessing the oscillating boundary given by the Legendre 
polynomials of various kinds [jj which lead to both reg- 
ular and chaotic motion pj. Similar investigations for 
the gravitational billiards were lacking. Only recently 
they have been studied - experimentally - under the as- 
sumption that the energy loss during the collisions is to 
be compensated, on the average, by the motion of the 



container [9j. The authors of Ref.[9j constructed three 
aluminum containers with parabolic, wedge, and hyper- 
bolic shape which exercised the horizontal oscillations. 
Inside a ball of steel was scattered from the boundaries 
and a camera registered the position of each bounce and 
the ball velocity. The results clearly indicate the regu- 
lar motion for the parabola and the chaotic one for the 
wedge; they also suggest some sort of regularity for the 
hyperbola at a small driving frequency. 

In this Letter we present a theoretical study of the 
inelastic gravitational billiards with a time-dependent 
driving. To the best of our knowledge, it is the first 
theoretical approach to these - very realistic - systems. 
The billiard shapes and parameters of the dynamical sys- 
tem have been so chosen to enable us a direct com- 
parison with the experimental data Then we as- 
sume the following boundaries: f(x) — a(x — x) 2 + c 
(the parabola), f(x) — b\x — x\ + c (the wedge), and 
f(x) = */ at(l + P(x — x) 2 ) — S (the hyperbola), where 
a = 0.26 cm" 1 , b = 1.85, c = 0.63 cm, a = 40.3 cm 2 , 
(3 = 0.08 cm~ 2 , and 6 — 4.45 cm. The containers oscil- 
late horizontally: x(t) = — Asimut, where A is the am- 
plitude and to — 2irf is the oscillation frequency. Inside 
the container the particle is subjected to the constant 
gravitational acceleration g. Collisions with the bound- 
aries result in the energy loss, quantified by the restitu- 
tion coefficient r 6 [0,1] which is defined as a ratio of 
the absolute values of the velocity after and before the 
collision. The case r = 1 corresponds to the elastic col- 
lision. It is difficult to decide a priori which value for r 
should be assumed. An experiment with steel particles 
bounced on a steel block gives r = 0.93 0. However, 
taking into account effects connected with the sharing 
of energy between rotation and translation during the 
collision reduces this coefficient substantially and the ef- 
fective r appears of about 0.7. Moreover, r can depend 
on the velocity and on the scattering angle [ljj. The 
authors of Ref.0 suggest r = 0.9. In the following, we 
will try to draw some conclusions about the restitution 
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coefficient from the comparison of our results with the 
experimental data of Ref.0- 

Let us assume that the particle hits the boundary at 
the time to with the velocity vfi , determined in respect 
to the frame connected with the billiard, and the collision 
point is (xq ,i/q — /(xq)). The transformation of parti- 
cle velocities at this point, — > vf , has the following 
form 

v? = r(v£ - 2u(v^ • u)), (1) 

where the components of the versor normal to the 
boundary, u, depend on Xq and are given by: u x — 

-f'{ x o)/h{ x o) and u v = l/H x o) witn H x o) = 
^Jl + f' 2 {x%). The particle, after being reflected from 
the boundary, moves along the ballistic trajectory: 

x c {t) = Xq + v%i(t - t ) + Asinujt - Asmujt 

y C (t)=Vo +v^ 1 (t-t )-g(t-t ) 2 /2, 

for t > to, up to the next section of this curve with 
the boundary: (x c (h) , y c (ti)) = (xf,j/f). The sub- 
sequent applying of Eqs.JTJ and @ produces a set of 
collision events which take place at times t n . The time 
evolution can be characterized by the vector X n = 
( x n >Vn i"i ' m v y n) m ^ ne four-dimensional phase space, 
where the velocities are taken just before the consecutive 
bounces. Therefore, we can restrict the dynamics to the 
billiard boundary and describe it in terms of the following 
time-dependent mapping: 

X„ +1 =P„(X„). (3) 

The above expressions have been formulated in the bil- 
liard coordinates because then the simple velocity trans- 
formation rule holds. The transformation to the 
laboratory frame is straightforward: x — x c — Asm cut 
and v x = v x — Auicosuit; y-components remain the 
same. The equation for the boundary f{x) reduces the 
4-dimensional phase space, in which V n is defined, to the 
3-dimensional manifold. 

For various shapes of the billiard and various driving 
forces, the mapping V n can represent either a regular 
cycle or a strange attractor. The static, elastic billiard 
with the parabolic shape is always regular and possesses 
two stable orbits: the horizontal orbit, connected with 
the symmetric bouncing between the left and the right 
part of the boundary at the same height, and the vertical 
one which involves the top of the parabola I n our 
case the limit cycle corresponds to the fixed point of the 
mapping Vl and can be obtained analytically by solv- 
ing the equation X* = "P 2 (X*). The detailed equations 
are complicated and will not be presented here. As a re- 
sult, we yield the stable horizontal orbit which moves to 
and fro together with the container. The vertical orbit 
of the elastic billiard shrinks now to a single point. The 
time interval between consecutive bounces of the horizon- 
tal orbit (2-point attractor on the boundary) r* = tt/uj, 



moreover v y = Trg/2uj; these quantities depend neither on 
the restitution coefficient r nor on the billiard shape. For 
the parabolic shape the agreement of t* with the data 
is very good. For r — 0.9 we get the height of the 
orbit y — 158 cm which exceeds by far the experimental 
value y « 7.5 cm. The latter value can be obtained if we 
assume r — 0.43. Therefore, the restitution coefficient 
seems to be well established by the experiment at the 
value r = 0.43 and we can try to apply it in numerical 
calculations for the other shapes. However, all trajecto- 
ries calculated in this way, both for the wedge and for the 
hyperbola, do collapse to the bottom of the billiard in- 
stantly. Apparently, the assumption that the restitution 
coefficient can be approximated by a constant cannot be 
maintained in the present problem. 

In the following we assume that the damping factor 
depends linearly on the particle velocity and then the 
restitution coefficient is of the form: 

r = 1 _| v C| /Vcr (4) 

if r > 0.01; otherwise (for very large velocities) r = 0.01. 
The parameter v cr can be determined by comparison 
with the data. The parabola case is especially useful for 
that purpose because it must be characterized by a sim- 
ple cycle. The data seem to confirm the presence of a 
cycle though it appears in a strongly diffused form which 
may result both from the low resolution in the experiment 
(e.g. due to the surface roughness) and from the rota- 
tional degree of freedom, not taken into account in the 
calculations. We get the position of the 2-point attractor 
in agreement with the data for v cr = 390 cm/s. We apply 
this value in the calculations for the other shapes. 

Fig.l presents some properties of the map V n - the 
height y and the time interval between consecutive 
bounces r„ = t n — f„_i - for the wedge and the hy- 
perbola which is driven by three different frequencies. 
Each figure represents a single trajectory, evolved up to 
t = 2 • 10 3 s which corresponds to about 2 • 10 4 colli- 
sions with the boundary. The plot of y for the wedge 
is strongly nonuniform and indicates a high degree of 
chaoticity. The most of the points is concentrated just 
above the line y n +i = Un and the fractal structure is 
hardly visible. The figure can be directly compared with 
the data (see Fig. 3 in Ref.0); the similarity is striking 
though the calculated quantities are extended to slightly 
larger values than the data. The nonuniformity is clearly 
visible also in the plot of the time intervals and the re- 
gion close to the point (0.09, 0.09) is distinguished. The 
experimental data exhibit the same pattern. The hy- 
perbola involves the wedge as well as the parabola as 
its limiting shapes and one can expect both shapes in- 
fluence the results of the dynamical calculations for the 
hyperbola. For the low driving frequency / = 4.5 Hz 
(Fig.lc, d) the particle abides close to the bottom of the 
billiard which can be well approximated by the parabola. 
That results in apparently regular pattern. However, a 
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FIG. 1: The coordinate y of a collision point vs. y for the 
previous collision (left column) and the time of flight between 
two subsequent collisions (right column) . The following cases 
are presented: the wedge with / = 6.6 Hz ((a) and (b)) and 
the hyperbola with: / = 4.5 Hz ((c) and (d)), / = 5.8 Hz ((e) 
and (f)), and / = 8 Hz ((g) and (h)). The units are: [cm] for 
y and [s] for the time. 



stochastic ingredient, similar as in Fig. la which corre- 
sponds to the wedge, is also visible in Fig.lc. In the ex- 
perimental data the chaotic component connected with 
the wedge seems to be absent completely. Instead, long- 
time tails in the plot of time intervals, possessing the in- 
termittent structure, are observed, as well as very small 
values of the height y. Then the experimental results 
suggest a more regular motion than the calculations pre- 
dict. For the larger frequency / = 5.8 Hz (Fig.le, f) 
the chaotic behaviour is overwhelming. The bands typi- 
cal for the strange attractors are visible; it is not the case 
for the experimental data but such subtle structures may 
be smeared due to the low resolution. The bands vanish 
completely if we make the frequency still larger. Fig.l 
g, h shows that it happens for / = 8 Hz (the case not 
studied experimentally) and the picture is similar to that 
for the wedge. However, it is not generally true that the 
degree of chaoticity rises with the frequency: for / > 8 
Hz the motion becomes regular again. 

Some additional information about the phase space 
structure can be obtained by plotting the position y ver- 
sus the tangential velocity after the collision, v t , at colli- 
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FIG. 2: The normalized coordinate y of a collision point vs. 
the normalized tangential velocity v t at the same point for 
the following cases: the wedge with / = 6.6 Hz (a) and the 
hyperbola with / = 4.5 Hz (b), with / = 5.8 Hz (c), and with 
/ = 8Hz(d). 



sion points, normalized by y max and w™ aa: , respectively. 
The quantities y max and v™ ax mean the largest possible 
values of y and Vt at each collision, obtained under the 
assumption that the particle is either at rest or bounces 
almost horizontally at the bottom of the billiard. There- 
fore, the region around v t — corresponds to the head- 
on collisions which may result in instability of the mo- 
tion and the onset of the chaotic behaviour. For the 
parabolic shape we have obtained 2-point limiting cy- 
cle, corresponding to very small values of the tangential 
velocity: \v t /v™ ax ,y/y max ) = (±0.044,0.39). The ex- 
perimental result presents itself as a narrow band, due 
to the noise in the system. Results of the calculations 
for the other shapes are presented in Fig. 2. The case of 
the hyperbola for the largest frequency / = 8 Hz shows 
the greatest disorder (Fig. 2d) whereas the trajectory for 
the wedge (Fig. 2a) is able to fill the entire region which 
is allowed by the energy partition rule |9j . The figure for 
the hyperbola at the intermediate frequency / = 5.8 Hz 
(Fig. 2c) indicates, in turn, a pronounced fractal struc- 
ture. On the other hand, the pattern for / = 4.5 Hz 
is predominantly regular (Fig. 2b). However, a vertical 
strip at small \vt \ bears apparent signs of chaos. Indeed, 
a magnification of this region reveals the fractal struc- 
ture. Also the experiment |9j distinguishes the region of 
small tangential velocities for this case but the lack of any 
fine structure in the data prevents detailed comparisons. 

A picture which has emerged so far is not complete yet 
because we do not have any quantitative knowledge - like 
the degree of instability and the attractor dimensionality 
- about the chaotic motion visible in the figures. Then 
we calculate the Liapunov exponents for the mapping V n 
which characterize the evolution of dynamical systems in 
the tangent space. Since the billiard is not smooth and 
the direct linearization of the equations of motion is not 
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TABLE I: The Liapunov exponents for the mapping Vn and 
the Kaplan- Yorke dimension for the parabola (p), the wedge 
(w), and the hyperbola with frequencies: f=4.5 Hz (hi), f=5.8 
Hz (h2), and f=8 Hz (h3). 



Ai A2 A3 A4 Dky 



p 


-0.408 


-0.579 


-0.579 


-2.17 


0.00 


w 


0.447 


-0.288 


-1.17 


-1.98 


2.14 


hi 


0.208 


-0.579 


-1.05 


-1.42 


1.36 


h2 


0.367 


-0.394 


-1.21 


-1.78 


1.93 


h3 


0.611 


-0.0637 


-1.05 


-2.08 


2.52 



possible, we resort to an approximate method [12,113 by 
utilizing the fact that the distance between two close tra- 
jectories is governed by the linearized equations. Then we 
perform the time evolution of two trajectories, whose ini- 
tial conditions differ by Sr, by the time interval St which 
must be small enough to keep the trajectories close. In 
the next step we renormalize the relative distance to Sr 
and continue the procedure for a long time. In order to 
get the entire spectrum of four Liapunov exponents Ai we 
need to evolve five trajectories and utilize the Gramm- 
Schmidt orthogonalization scheme at each renormaliza- 
tion step to sort manifolds connected with subsequent 
unstable and stable directions 0]. Finally, we have to 
multiply the exponents by the mean time between suc- 
cessive bounces, different for each case. The results are 
summarized in Table I. All presented cases are charac- 
terized by one positive exponent, except the parabola 
for which all exponents are negative. For the hyperbola, 
the instability rises with the frequency and even the tra- 
jectory for the smallest frequency is chaotic despite the 
apparent regularity visible in the figures. 

Having all Liapunov exponents calculated, we can 
determine the Kaplan- Yorke dimension which in many 
cases may be identified with the Hausdorff dimension 
(the Kaplan- Yorke conjecture ^5|). The definition is the 
following 

Dky = i + ^ i , (5) 

where j is the largest integer such that Ai + A2 + • • ■ + \j > 
0. The results, presented in Table I, indicate that the 
embedding dimension of the attractor for the cases w 
and h3 equals 3, i.e. it is equal to the entire available 
manifold. This conclusion explains the nature of chaotic, 
but also nonuniform and diffused, pattern observed in the 
figures for those cases, in respect both to the theoretical 
predictions and to the experimental data: for so high 
dimensionality of the attractor its structure cannot be 
clearly visible in the plane. For the wedge that diffused, 
nonfractal pattern is restricted to the area just above the 
line y n +i = Vn but in the case h3 it extends to the almost 



whole picture. On the other hand, the fractal structure 
is apparent for the case h2. 

We have demonstrated that the inelastic gravitational 
billiard with the time-dependent driving, possessing the 
shape of the wedge and the hyperbola, constitutes the 
strange attractor, whereas the parabolic shape is charac- 
terized by the regular motion. We have proposed a simple 
dependence of the restitution coefficient on the velocity. 
Such dependence appears to be essential to get results 
consistent with the experimental data. Some complicated 
patterns, revealed by the experiment, can be explained 
by the attractor dimensionality and its dependence on 
the billiard shape and parameters. Generally, the model 
predictions agree quite well with the data. However, the 
values of the height y are slightly too large, which re- 
sults in the more irregular motion than the experiment 
shows, for the hiperbolic shape with small driving fre- 
quency. This discrepancy may be a token of too weak 
damping and a suggestion that the simple formula (0J 
should be refined at small velocities. For that purpose 
some experimental effort is necessary: the resolution of 
the data should be improved and, first of all, the resti- 
tution coefficient, which is the essential quantity in the 
model, should be determined in the wide range of the 
ball velocity. 
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